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LIST OF SYMBOLS 
Definition 

a quantity subscripted by T pertains to the entire system 

a quantity subscripted by P pertains to the instrument (body P) 

a quantity subscripted by G pertains to the gondola (body G) 

a dot over a quantity denotes its time derivative 

an overbar denotes a vector quantity 

the acceleration due to gravity 

mass 

weight 

unit vectors — the subscript a denotes the coordinate frame and the 
subscript b denotes the axis direction, b = 1,2,3 

coarse elevation angle traversed by the coarse elevation gimbal axis 

Euler angles - subscript a denotes the coordinate frame and subscript 
b denotes the axis about which the angle is measured, b = 1,2,3 

magnitude of the Euler angle rates (first time derivatives of 6 ab ) 

acceleration of the balloon attach -point 

position vector of the center-of-mass of body P 

position vector of the center-of-mass of body G 

position vector of the center-of-mass of the system 

force in the cable f^ = f^ i^ 

vector representing the cable = Iq i^ 

force on body G due to gravity (weight of G), f G ^ = " Wq i Lgl 

force on body P due to gravity, f p ^ = - W p i LS1 

vector from_ origin of C frame to center-of-mass of body G, 

P G = " P G l Gl 

vector from_ origin of C frame to center-of-mass of body P, 
p P = ‘ p P : G1 

vector from origin of C frame to center-of-mass of the system, 

P CM = " P CM i Gl 
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p G-P 
IF 


ext 


ZM 


CM( ) 


denotes the scalar quantity p G - Pp 

sum of external forces acting on the system 

sum of moments about the center- of- mass of a body (the character in 
the subscript parenthesis denotes the body) 


H CM( ) 


ZF 


(P) 


“P 

“G 

T 

P/G 


T 

A G/C 

k g/c 

f 

P/G 


angular momentum about the center-of-mass of a body 
sum of the forces acting on body P 

angular rate of body P, uj p = u pi i pi + oj p2 f p2 + u) p3 f p3 
angular rate of body G, = a> G1 i G1 + u Q2 i Q2 + co G3 f G3 
torque on body P due to body G, (note T p ^ G = - T G y p ), 

T P/G = T P/G(P1) *P1 + T P/G(P2) 1 P2 + T P/G(P3) *P3 

torque on gondola due to cable torsion, T G ^ C = - K G / G 0 G1 i G1 

torsional spring constant of the cable 

force on body P due to body G, (note f p ^ G = ~ ^ G /p)» 

f P/G = f P/G(LSl) 1 (LSI) + ^P/G(LS2) X LS2 + f P/G(LS3) ^83 




X 

tJ 

9j, i=l, 2, 3 
K Ri , i=l,2,3 


K p ., i=l,2,3 


output torque of the azimuth torque motor 

output torque of the fine elevation torque motor 

output torque of the roll torque motor 

resistance torque due to azimuth bearing friction 

resistance torque due to fine elevation axis flex pivot 

resistance torque due to roll axis flex pivot 

azimuth torque less bearing friction 

elevation torque less flex pivot stiffness effects 

roll torque less flex pivot stiffness effects 

denotes vector of state variables 

denotes vector of input variables 

pointing angles — Euler angles between the P and S reference frames 
rate gain in the control law 
proportional gain in the control law 


vi 



K n , i=l, 2,3 


G 


integral control gain in the control law 
inertia tensor for the P body 
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a ph 

*P12 

J P13 

*p 

= 

J P21 

X P22 

X P23 



T P31 

J P32 

! P33 

inertia tensor for the G 



l Gll 

J G12 

X G13 

*G 

= 

I G21 

! G22 

X G23 



*G31 
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! G33 


vii 



SUSPENSION 
SHROUD LINES 













A universal joint is located between the gondola and the balloon suspension cable 
to reduce bending moments between the gondola and the cable. 

Using the gondola and gimbal system configuration outlined above, a math model 
of the system can now be developed. 


MODELING THE SYSTEM 


The system is modeled as two rigid bodies connected together by a gimbal sys- 
tem. Figure 2 shows the simplified model and the coordinate systems used to derive 
the equations of motion. The instrument and its counter- weight form body P, which 
is connected through the gimbal system to the gondola structure (body G). The 
gondola is suspended from the balloon- attach-point by a cable. 
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Figure 2 




The LS, S, and B coordinate systems all have their origins fixed at the balloon- 
attach-point. The LSI axis points away from the Earth along the local vertical. The 
LS3 axis is directed so that the Sun lies in the LSI to LS3 plane. The SI axis 
points to the Sun and the S2 axis is aligned with LS2. 

The B coordinate system is fixed in the balloon and rotates with it. Only 
balloon rotation about LSI is considered, so B1 remains aligned with LSI as the 
balloon rotates. 

The C coordinate frame has its origin fixed at the universal joint where the 
cable and gondola connect. The Cl axis remains aligned with the cable. The origin 
of the G coordinate frame is at the center- of- mass of the gondola and the G1 axis 
points along the vertical pylon of the gondola as shown. The origin of the P coor- 
dinate system is fixed at the center- of- mass of the P body. The Pi axis is aligned 
with the line-of- sight of the instrument. 

Any friction in the balloon-cable and cable-gondola connections is neglected, as 
is any internal cable friction. The dynamics of the balloon are limited to a transla- 
tional acceleration, a„, applied to the balloon-attach-point, and an angular rotation 
about the B 1 axis . 

Several other general assumptions apply to the model. Any aerodynamic dis- 
turbances on bodies P and G are neglected. The cable has a torsional stiffness, but 
is assumed to be massless and extensionless. The center-of-mass of each body lies 
along the G1 axis. Finally, it is assumed that the dynamics of the cable, body P, 
and body G do not affect the balloon motion. 


DERIVING THE EQUATIONS OF MOTION 


The orientation of the system at any time can be described by the Euler angles 
between the coordinate systems. Two angles are needed to express the transforma- 
tion from the B frame to the C frame, 0 c2 and 0 c3 . Three angles are required from 

the C frame to the G frame, © G1 , 0 Q2 , and 0 G3 * Finally, 0 pi , 0 p2 , and 0p 3 define 

the transformation from the G frame to the P frame. Eight equations of motion, con- 
taining these angles and time derivatives of these angles, are required to express the 
motion of the system relative to the B coordinate frame. 

Several more assumptions are made in order to simplify the equations. First, 
it is assumed that the Euler angles will be small, so that cos 0 % 1 and sin 0 % 0. 
Also, higher order terms in products of state variables and/or inputs are neglected. 
During the initial derivation of the equations, the azimuth bearing friction and the 
flex pivot stiffnesses are neglected. The friction and stiffnesses will be modeled and 
included in the equations later in the derivation. 

These equations can be derived in four basic steps. The first step is to solve 
the equation 


ZF 


ext 




+ a B ) 


( 1 ) 


3 • 



r^ is defined relative to the LS coordinate frame. Since this frame accelerates with 
respect to inertial space, the addition of the term § B is necessary. The external 
forces acting on the G body and the P body are the cable force, , the gondola 
weight, ?Q(gy» and the instrument weight, f p ^^ ; therefore. 


ZF ext ~ f C " f G(g) " f P(g) 


( 2 ) 


In terms of unit vectors, equation (2) can be written as 


ZF ext f C *01 “ (W P +W G ) ^Sl 


(3) 


The vector r^ can be expressed in terms of the physical dimensions z c and p CM as 


r CM - ' Z C 1 C1 " P CM 1 G1 


(4) 


From this expression r^ M can be obtained, 
point, a fi , can be expressed as 


The acceleration of the balloon-attach- 


a B ^(LSl) ^Sl + a B(LS2) 1 LS2 + ^(LSS) *LS3 


(5) 


Notice that the components of equations (3), (4), and (5) are expressed in different 
coordinate systems. Coordinate transformations are needed to express all of the com- 
ponents in a common coordinate frame. Transforming all of the quantities in the 
above equations to the LS frame, and substituting the results into equation (1), the 
following equations result. 


f C - W T + m T a B(LSl) 

~ m T (p CM +£ C* 8 C3 " m T P CM ®G3 


(6a) 


(W T + m T a B(Lgl) ) 6 C3 


* n T a B(LS2) 

(6b) 


m T (p CM +£ C ) 


d C2 + m T P CM 



(W T + m T a B ( LS1 )) 0 C2 


ra T ®B(LS3) * 


(6c) 
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Equation (6a) gives an expression for the force in the cable. Equations (6b) and 
(6c) are equations of motion containing Euler angles and second time- derivatives of 
the Euler angles. 

The next step is to sum the moments about the center-of-mass of the P body 
and equate this to the time derivative of the angular momentum, 


ZiV1 CM(P) H CM(P) 


(7) 


The time derivative of the angular momentum is given by the expression 

H CM(P) = *p “p + w p x ^ 

where Ip is the inertia tensor for the P body and. u p is its angular velocity. The 

term co p x (I p a> p ) consists entirely of second-order quantities and can be neglected 

since higher-order terms are assumed to be negligible. The sum of the moments 
acting on body P is the torque applied to the instrument by the gimbal torque motors. 


zM CM(P) t p/g 


(9) 


To express T p ^ G in terms of the gimbal axis torques, the transformation from the G 

coordinate system to the P coordinate system must be considered. Appendix A shows 

• • 

how T p ^ G is expressed in terms of the gimbal torques. Substituting = I p w p 

and the expression for T p/ G from Appendix A into equation (7), the following three 
equations of motion are obtained. 


! pll “pi + *pl2 “p2 + ! pl3 “P3 = T A 008 0 CE " T E* Sin 0 CE 

*p21 “pi + X p22 “p2 + *p23 “p3 = T e' 

*p31 “pi + *p32 “p2 + J p33 “p3 = T A* sin 9 CE + T r' 008 0 CE 


(10a) 

(10b) 

(10c) 


The above equations do not contain any Euler angles or their time derivatives expli- 
citly, but they could be obtained through the relationship between ui p and the Euler 

angle rates. The equations are retained in the above form because the rates u^, 

Upg, and u p g will later be used as state variables. They are quantities that are 

measured by the rate gyros. 
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The third step in deriving the equations of motion is to solve the vector 
equation 


y M = H 

ilv CM(G) n CM(G) 

As for the P body, 

Z ^CM(G) = " ^G x 


( 11 ) 

= IqCOq . The moments acting on the G body are given by 

(P G -Pp) x f G /p + T q/p + f Q/c . (12) 


The torque T G ^ p is the torque applied to the gondola by the gimbal torquers. It is 

equal to the negative of T p ^ G derived in Appendix A. ^q/q is the torque applied 

to body G due to the torsion in the cable. t q/ G is derived in Appendix B. The 

force f G y p is the interconnection force between the G body and the P body. When 

summing moments on body P, this force did not cause a moment because the center-of- 
mass of the P body is at the center of the gimbals, which is the point of application 

of the force f G / p * This force does, however, cause a moment about the center-of- 

mass of the gondola. By transforming each element of equation (12) into the G coor- 
dinate frame and substituting these results into equation (11), three more equations 
of motion can be obtained. However, one element in equation (12) is still unknown. 

The force f G / p must be determined. By solving a fourth basic equation, the force 

fp/G , which is the negative of ? G / p » can be found. Solving the equation 


zF (P) ~ m P ^ r P + a B^ 


(13) 


will yield the force The sum of the forces on body P can be expressed as 


EF (P) V/G + ? P(g) 


(14) 


where ^p(g) is simply the weight of body P. As with equation (1), all the elements 

in equation (13) are transformed to the LS frame. This results in three equations 

that give the three components of the force f p / G > which is the negative of f G / p 

needed for equation (12). These expressions are transformed to the G coordinate 
frame where they are substituted into equation (12). The resulting equations of 
motion from equation (11) are 


X G11 “G1 + ! G12 6 G2 + ^lS “G3 = " T a' " K G/C 0 G1 


(15a) 
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l G21 “G1 + X G22 W G2 + *G23 W G3 = ' T e’ “ P G m T (g + a B(LSl) ) 6 G2 

+ ( Pq~ p p) [( 9 G2 +6 C2^ ^ W p +m p a B(LSl)^ 

+ m p (SL C e c2 + a B(LS3) + m P p P ^G2 + ^C2^ 

(15b) 


*G31 U G1 + *G32 W G2 + *G33 W G3 = " T R* " P G m T (g + a B(LSl) ) 0 G3 

( p G -p p) ^'^ 6 G3 +e C3^ ^ W p +m p a B(LSl)^ 

+ m p ( a B(Lg 2 ) " ®C3^ ' m P p P *®G3 + ®C3^ 

(15c) 

The above equations contain the components of the angular rate vector oj^. 

These components can be eliminated from the equations by expressing them in terms 
of the Euler angle rates. Appendix C shows how this substitution is derived. 


DEFINING THE STATE VARIABLES 


Most of the equations of motion that were derived are second-order differential 
equations. These equations can be transformed to first-order equations by defining 
eight additional state variables. Let fi denote the time rate of change of an Euler 
angle , or 


n = e 


(16) 


Then the time derivative of ft will be 


n = 0 


(17) 


and the substitution can be used to transform the equations of motion to first-order 
form. Equation (16) will give additional state equations needed to find the time- 
histories of the Euler angles 0 G2> 0 G3> © G1 , 0 G2> and 0 G3 . 
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Since was not eliminated from equations (10a), (10b), and (10c), equation 
(16) cannot be used to express the time histories of 0 p ^, 0 p ,, , and 0 p . The relation- 
ship between tu p and the time derivatives of these three Euler angles must be derived. 

The details of developing this relationship are explained in Reference 3. The final 
three state equations that result from this relationship are 


Qpi COS 0 CE + U)p 3 sin 0 CE - n G1 ft B1 


(18) 


0p 2 w p2 " fi C2 


(19) 


e p 3 --wpi sin 0 CE + o) p3 cos 0 CE - fi G3 - fi c3 


( 20 ) 


The other 13 state equations previously derived, with all substitutions and sim- 
plifications included, are given below. 

m T (p CM + £ C ) ^C3 “ m T P CM ^3 = (W T + m T ^(LSl)^ 6 C3 " m T a B(LS2) 

( 21 ) 

m T (£ C + p ClVp ^C2 + m T P CM ^G2 = " (W T + m T a B(LSl) ) 0 C2 “ m T a B(LS3) 

( 22 ) 


r Pll “PI + J P12 W P2 + ! P13 “P3 T a' COS e CE " T r’ sin 0 CE 


(23) 


*P21 W P1 + *P22 a) P2 + ^23 “P3 T E ? 


(24) 


J P31 ^Pl + *P32 “P2 + r P33 W P3 T A* sin 0 CE + T r’ COS 0 CE 


(25) 


^11 fi Gl + I G12 Q G2 + r G12 n C2 + *013 ^G3 + I G13 ^C3 _ ' T a' K G/C 0 G1 


!gii n B1 (26) 
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®C2 “ n C2 (29) 

®C3 = ^C3 (30) 

®Gl = fi Gl (31) 

®G2 = fi G2 (32) 

®G3 = fi G3 ' (33) 


The eight second-order differential equations have now been replaced by 16 first-order 
differential equations. Using matrix notation, these equations are in the form 


MX = KX + CU . (34) 
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This can be easily transformed to the common state variable form. 


X = A X + B U 


( 35 ) 


by inverting the M matrix. This inversion is easily performed in the computer simula- 
tion, so it is not necessary to perform it analytically. 

The state vector X is defined as 


fi C3. 

°Q1 

fi G2 

^G3 

“PI 

u) p2 




i '» 


P3 

0 C2 

0 C3 

0 G1 

0 G2 

0 G3 

0pi 

®P2 

^3 


> 




J 


(36) 


It contains the eight Euler angles that were originally sought to define the position of 
the system. The Euler rates of the C and G coordinate transformations, and the body 
rates of the instrument, are also state variables. 

The input vector U is defined as 
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(37) 


U = 


' t a' ^ 

T ' 

E 


l R 

S a B(LS2) 


a 


B(LS3) 




B 1 


B1 




which contains the gimbal torques and the balloon motions. 


Notice that the A matrix in equation 
the balloon has a nonconstant acceleration, 
is present in the A matrix. 


(35) will be a function of time whenever 
This is because the component 


Until now, the azimuth bearing friction and the stiffness of the flex pivots have 
been neglected. Modeling these parameters will now be considered. 


MODELING THE GIMBAL BEARINGS 


In the present effort, a very simple model is used for the gimbal bearings. 
This model is sufficient for the preliminary control system design, but a more accu- 
rate model may be desired if a more precise analysis is needed. 

The azimuth bearing is modeled using a Dahl friction model, which is explained 
in detail in Reference 2. Using small angle assumptions and neglecting higher-order 
terms, the following equation for the resistance torque of the bearing results. 


T 


fA 




(38) 


The parameters y and T^ q are constants that depend on the physical properties 

of the bearing. Using values given in Reference 2 as guidelines, values are chosen 
for y and T^ Q that will yield reasonable values for T^ relative to the azimuth gimbal^ 

torque, T A « Again, for a more precise analysis, these values would have to be 

determined for the particular azimuth bearing. 

Equation (38) can be integrated directly to give the resistance torque 

T fA = J A 0 P1 (39) 
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where is a constant, 


J 


A 



The flex pivot bearings act as torsional springs and are, therefore, modeled as 
linear torsional springs with a resistance torque given by 


T 

res 



e 


(40) 


where J is the spring constant. Values for J are available from Reference 1. The 
resistance torques of the two flex pivots are expressed as 


T fE J E 0 P2 

(41) 

T fR = J R 0 P3 

(42) 


for the fine elevation and roll axes, respectively. To include the bearing models into 
equations (18) through (33), the following substitutions are made in those equations 


V = t a + J A 6 P1 


(43) 


V T E + J E e P2 


(44) 


V " T R + J R ®P3 


(45) 


where T^, Tg, and T^ now represent the actual torque applied by the torque motors 
on the gimbal system. 

This completes the development of the equations of motion. The equations^ can 
now be used to simulate the motion of the system for the various inputs in the U 
vector. The next step is to develop a control law to point the instrument at the Sun. 
This control law will determine the gimbal axis torques required to point the instru- 
ment. 
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DEVELOPING THE CONTROL LAW 


Since the objective is to point the instrument at the Sun, a relationship between 
the S reference frame and the P reference frame must be developed. If it is assumed 
that azimuth and coarse elevation axes motions have pointed the instrument relatively 
close to the Sun, the small angle assumption can be used when defining a 
transformation between the S and the P reference frames. Using this assumption and 
neglecting the second-order terms, an Euler transformation can be developed. 


r 






> 

1 pi 


1 

9 3 

1 « 
CD 
1 


C _ 

1 S1 

1 P2 

> = 

‘ e 3 

1 

9 1 

< 

1 S2 

1 P3_ < 


_ 0 2 

- 0 1 

1 


1 S3 

L J 


(46) 


Another transformation matrix can be developed between the S and P frames by 
making successive Euler transformations from the S to LS, LS to B, B to C, C to G, 
and G to P coordinate frames. This transformation matrix will be a function of the 
Euler angles of the system that result from the equations of motion. By equating this 
matrix to that of equation (46), the angles 0 1 , 6 2 , and e 3 can be determined in terms 

of the other Euler angles. The following expressions are found for these angles. 


9 1 = (0 B1 + 6 G1 + 0 Pp COS 0 CE “ (0 C3 + 0 G3 + ^3* Sm 0 CE (47) 

0 2 = 0 C2 + 0 G2 + 0 P2 (48) 

0 3 = ^ 0 C3 + 0 G3 + 0 P3^ COS 0 CE + ^ 0 B1 + 0 G1 + 9 Pp sin 0 CE * * 49 ^ 


Angles 0p ©2» and 0^ are the Euler angles between the P and S reference 

frames. When these angles are zero, the two frames will be essentially colinear and 
the instrument will be essentially pointing at the Sun, because the distance to the 
Sun is much greater than the distance between the origins of the P and S reference 
frames. Now that the pointing angles 0^, ©g, and 0g have been determined in terms 

of the Euler angles of the system, a control law can be developed which will maintain 
these angles close to zero. 

The control law for each axis is defined to be 


^(Pi)com 


K Ri [K Pi (0 (i)com 


- V - 


a pi 


■Vii 


(50) 
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The command angle 0 ^) com is included so that the instrument can be pointed 

at a source slightly different from the center of the Sun. Figure 3 shows a single 
axis block diagram of the control law along with a simplified model for the plant 
dynamics . The intergral control is required to remove any steady- state offset in the 
response. The rate commands w„ T and must be limited due to the limited control 

torque T 


cl 


cP 


(Pi)com* 


These limiters are shown in the block diagram. The limiter on 


a j must be such that the integrator stops integrating while u)^ is at its limit. 


The commanded torques obtained from equation (50) are the commanded torques 
to be applied to the instrument along the axes of the P coordinate system (body axes). 
These commanded torques must be transformed as shown in Appendix A to determine 
the commanded gimbal torques T (A)com , T (E)com , and T ( „ )00(n . 

In the present analysis, the gains K^, Kp., and Kp. were chosen to yield a 

phase margin of 65 deg and a gain crossover frequency of 1 Hz for the system in 
Figure 3 when linearized. 



Figure 3 


SIMULATION RESULTS 


A computer simulation of the system was written using equations (18) through 
(33) and the control scheme developed in the previous section. The sampling rate of 
the digital flight computer, 10 samples per second, was included. The sensors and 
torque motors were assumed to be ideal, except for a limit on the maximum torque 
available from the torquei motors. A listing of the simulation code is included in 
Appendix D. 
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Simulation results are shown in Figures 4 through 8. These results were 
generated with zero initial conditions for all of the state variables. This means the 
instrument is initially pointing at the Sun. A step input for the command angles 6^, 

0 2> and 0g of 0.003 radians is introduced at time t = 0. A limit of 0.0006 rad/sec in 

each axis was placed on the rate command in the integral control path. A limit of 
this same value was also placed on the rate command in each axis in the proportional 
control path. These limits resulted in a desirable response. Figures 4, 5, and 6 
show the gimbal axis torques and the response of the angles 0^, 0 2> and 0^ to the 

step input. The elevation and roll torques saturate at times. Figure 7 shows the 
components of the angular rate of the instrument. 

Figure 8 shows the components of the pendulous motion of the gondola rela- 
tive to the C coordinate frame. This behavior is typical of the cable motion relative 
to the balloon also. 


CONCLUSIONS 


The balloon-pointing- system was modeled as a two rigid-body system. Equa- 
tions describing the motion of the system were derived and a control law was adopted 
to point the science instrument at the Sun. A simulation program was developed and 
results from the simulation indicate that the instrument can be pointed with favorable 
response characteristics . 

Some simplifying assumptions were made during the derivation of these equa- 
tions of motion. The equations, however, are sufficient for a practical analysis of 
the control system. 
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APPENDIX A 


The commanded torques from the control law are commanded torques in the body 
axes of the instrument. That is, 


T P/G T P/G(P1) X P1 + T P/G(P2) X P2 + T P/G(P3) X P3 * 


(A-l) 


The commanded gimbal torques, however, are commanded torques along the Euler 
axes of the transformation between the G and P frames. Four rotations are used to 
define this transformation: azimuth, fine elevation, fine roll, and coarse elevation. 

The transformation is expressed mathematically as 


/■ 


■ — 

X P1 


cos 0 CE 

) X P2 

v = 

0 



M., 

CO 

\ 


_ sin 6 CE 


where small angles are 


0 

1 

0 


-sm 0, 
0 

COS 0 


CE 


CE 


0 . 


P3 

0 


0p3 0 


1 0 
0 1 


1 

0 

^2 


assumed except in the case of the 


■ —i 


- 


f . ' 

0 -0p 2 


10 0 


X G1 

1 0 


0 1 Op! 


X G2 7 

0 1 


0 -ep! 1^ 


X G3 


(A- 2) 


coarse elevation angle 0 cg . 


The commanded torque on body P can be expressed in terms of the commanded 
gimbal torques for the azimuth, fine elevation, and fine roll axes as. 


P/G 


= T 


G1 


+ T, 


P2 


+ T 


R P3 


(A-3) 


The coarse elevation gimbal is assumed to be mechanically locked in this develop- 
ment. The primes are put on T^, T.g, and T R here because the effects of the gimbal 

bearings have not yet been iricluded in the model. This notation will avoid confusion 
later on. The unit vector ip 2 * is the unit vector along the new 2-axis that results 

after the rotation of 0p. about the G^ axis in equation (A- 2). The unit vector ipg” 

is the new 3- axis that results from the rotation through angle 0p 2 * Using equation 

(A-2), the relationships between i^, i pg f » ipj", and the unit vectors of the P 

reference frame can_be determined.. Equation (A-3) can then be expressed in terms 
of the unit vectors ip,, i p 2 , and f, g. Equating this result to equation (A-l) gives 

the relationship between the control torques and the gimbal. axis torques. These 
relationships are, in matrix form, 
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cos 0 CE 

_0 P3 

sin 0 CE 



' ' 

T P/G(P1) 

0 P3 cos 0 CE 

1 

0p 3 sin 

0 CE 


T P/G(P2) 

-(sin 0 CE + 0 p2 cos 0 CE ) 

0 

(cos 0 CE - 

0p 2 sin e CE ) 


T P/G(P3) 


(A- 4) 
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APPENDIX B 


The torsion in the cable can be expressed as 


t g/c " k g/c a 'ci 


(B-l) 


where a is the angle of twist in the cable. For the small angle assumption, the angle 
a will be about the same as the Euler angle 0^, 


a ^ 


G1 


(B-2) 


From the Euler transformation between the C and the G frames, i^ can be expressed 
in terms of the G frame unit vectors as, 


1 C1 1 G1 " 6 G3 1 G2 + 0 G2 *G3 


(B-3) 


assuming small angles for and Substituting equations (B-2) and (B-3) into 

equation (B-l) and neglecting second-order terms, the cable torsion becomes 


T G/C “ “ K G/C e Gl *G1 


(B-4) 
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APPENDIX C 


The vector is the angular rate of the gondola. Referring to Figure 2, this 


can be expressed as 

“a = S B1 *B1 + ®C2 l B2 + Q C3 1 C3 + 

+ ®G1 1 G1 * 

Let 

&B1 = n Bl 

e C2 = fi C2 
* 

0 C3 = fi C3 
®G1 = n Gl 
0 G2 = Q G2 
0 G3 = ^G3 


®G3 1 C3 + ®G2 (cos 9 G1 1 G2 ' sin ®G1 1 G3 ) 

(C-l) 

(C-2) 

(C-3) 

(C-4) 

(C-5) 

(C-6) 

(C-7) 


Substituting equations (C-2) through (C-7) into equation (C-l), expressing the 
unit vectors in equation (C-l) in terms of the unit vectors Oq^> i G2 > *G3^’ assum ^ n S 


sin 6=0 and cos 0 = 1 for all the 9’s in the resulting equation for (C-l), and neglect- 
ing the higher order terms in the 0's, the 0's and their product yields 


^G W G1 'Gl + “G2 1 G2 + W G3 1 G3 


(C-8) 


where 


W G1 fi Gl + fi Bl 


(C-9) 
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(C-10) 


“G2 fi G2 + n C2 


“G3 “ Q G3 + fi C3 


(C-ll) 


Taking the time derivative of equation (C-8) and neglecting higher order terms in 
a) G ., i = 1,2,3 and their product produces 


“G ~ “G1 1 G1 + “G2 1 G2 + W G3 ^3 


(C-12) 


where 


W G1 “ fi Gl + ft Bl 


(C- 13) 


“G2 “ a G2 + fi C2 


(C- 14) 


U G3 fi G3 + fi C3 * 


(C-15) 


Equations (C-13)S through (C-15) follow from equation (C-9) through (C-ll) 
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APPENDIX D. LISTING OF COMPUTER PROGRAM 


C++++++++++++4-+4-+-M- 4-++4-++'f4-++++++'f 


C PROGRAM TO SIMULATE P.E.'S PROPOSED GONDOLA CONCEPT AND + 
C GIMBAL SYSTEM FOR THE BALLOON POINTING SYSTEM + 
C THIS CODE CONTAINS THE SIMULATION OF THE EOM'S WITH THE ♦ 
C CONTROL LAW, INTEGRAL CONTROL, AND THE SAMPLE & HOLD OF THE + 
C DIGITAL COMPUTER INCLUDED. ALSO, LIMITS ARE IMPOSED ON THE + 
C TORQUE OF THE GIMBAL TORQUERS , THE RATE COMMAND DUE TO THE + 
C INTEGRAL CONTROL, AND THE RATE COMMAND DUE TO THE ANGLE ERROR. + 
C IMSL ROUTINES ARE USED FOR MATRIX INVERSION, MATRIX + 
C MULTIPLICATION, AND INTEGRATION OF THE STATE EQUATIONS. + 
C BILLY LIGHTSEY + 
C POINTING CONTROL SYSTEMS BRANCH/CONTROL SYSTEMS DIVISION + 
C APRIL 1988 + 


C+-++++++++++++++++++++++++++++++ ++++++++++++■►+•*-+++++++++++++♦++++++ 
IMPLICIT REAL(I.K) 

REAL M(16, 16) ,A(16, 16) ,B(16,7) ,X(16) ,U(7) ,Z(16, 16) ,MT,LC,MP,MG, 
ill ,12, 13, IM1 , IM2, IM3, UMAX, I2MAX, I3MAX 
DIMENSION PARAM(50) 

INTEGER IDO 
EXTERNAL FCN.IVPRK 

COMMON /SUB1/ A, B,2 , WT,MT , RG, RGP , WP , U, AB1 ,MP 

COMMON /SUB2/ AB10,AB11,AB20,AB21 , AB30 , AB31 , OMBO.OMB1 , OMGAB, TB10 
6TB 1 

COMMON /TIMES/ TAB 1 1 , TAB12 , TAB21 , TAB22 , TAB3 1 , TAB32 , TOMB1 , T0MB2 
NAMELIST /INPUT/ LC ,RP , WG.KGC.RG, WP, TCE, IP1 1 , IP12, IP13, 

&IP21 , IP22, IP23, IP31 , IP32, IP33 , IG11 , IG12, IG13, IG21 , IG22, IG23, IG31 
&IG32, IG33.AB10.AB11 ,AB20, AB21 ,AB30, AB31 ,OMBO, OMB 1 , OMGAB, TB10, 
&DELT, TI, TF, TABl 1 , TAB12 , TAB21 , WC, PM1 ,PM2, PM3, IPRT.X, 

LTAB22 , TAB31 , TAB32 , TOMB 1 , T0MB2 , TFA , TFE , TFR , TCI , TC2 , TC3 , TQMAX , 
&WCT1 ,WCT2,WCT3,WCI1 ,WCI2,WCI3 
READ (5, INPUT) 

MG-WG/32.2 
MP-WP/32.2 
MT-MG+MP 
WT-WG+WP 
CTCE-COS (TCE) 

STCE-SIN (TCE) 

RCN- (MP*RP+MG*RG) /MT 
RGP-RG-RP 

C CONVERT PHASE MARGINS TO RADIANS 
PM1-PM1*3. 14159266/180. 

PM2-PM2*3. 14159266/180. 

PM3-PM3*3. 14159266/180. 

C COMPUTE CORNER FREQUENCIES TO CORRESPOND TO THE PHASE MARGINS 
WC2-WC**2 
TNP-TAN(PMl) 

W1 1 -WC* (SQRT (TNP*TNP+ 1 ) -TNP) 

TNP-TAN (PM2) 

W12-WC* (SQRT (TNP*TNP+ 1) -TNP) 

TNP-TAN (PM3) 

W13-WC* (SQRT (TNP*TNP+1) -TNP) 

W21-WC2/W11 

W22-WC2/W12 

W23-WC2/W13 

C COMPUTE THE GAINS KI,KP,KR 

KIl-(WC2*SQRT(l+WC2/W21**2)) /SQRT (1+WC2/W1 1**2) 

KI2- (WC2*SQRT (1+WC2/W22**2) ) /SQRT (1+WC2/W12**2) 

KI3- (WC2*SQRT (1+WC2/W23**2) ) /SQRT (1+WC2/W13**2) 

KP1-KI1/W11 

KP2-KI2/W12 



KP3-KI3/H13 

KRI-W21*IP11 

KR2-W22*IP22 

KR3-W23*IP33 

C INITIALIZE M,A, AND B MATRICES TO ZERO 
DO 1 N-1,16 
DO 2 J-1,16 
M(N, J)-0.0 

2 A(N,J)-0.0 
DO 3 L-1,7 

3 B(N,L)-0.0 

1 CONTINUE 

C CREATE THE M MATRIX 

M ( 1 , 2) — MT* (RCM+LC) 

M(l,5)— MT*RCM 
M(2, 1)"MT* (LC+RCM) 

M(2,4)-MT*RCM 

M(3,6)-IP11 

M(3,7)-IP12 

M(3,8)-IP13 

M(4,6)-IP21 

M(4,7)-IP22 

M(4,8)-IP23 

M(5,6)“IP31 

M(5,7)-IP32 

M(5,8)-IP33 

M(6,3)-IG11 

M(6,4)-IG12 

M(6, 5) “IG13 

M(6,1)-IG12 

M(6,2)-IG13 

M(7,3)-IG21 

M (7 , 4) -IG22-RGP*MP*RP 

M(7,5)-IG23 

M(7, 1) ”IG22-RGP*MP* (LC+RP) 

M(7,2)“IG23 

M(8,3)-IG31 

M(8,4)-IG32 

M (8 , 5) -IG33-RGP*MP*RP 

M(8,1)-IG32 

M (8 , 2) -IG33-RGP*MP* (RP+LC) 

DO 4 N-9,16 

4 M(N,N)-1.0 

C INVERT THE M MATRIX 

CALL LINRG(16,M, 16, Z, 16) 

C CREATE NEEDED ELEMENTS OF K MATRIX. THE TIME-DEPENDENT ELEMENTS ARE 

C CALCULATED IN SUBROUTINE FCN 
K611— KGC 
K146-CTCE 
K148-STCE 
K166— STCE 
K168-CTCE 
K314— TFA*CTCE 
K316-TFR*STCE 
K415— TFE 
K514— TFA*STCE 
K516— TFR*CTCE 
K614-TFA 
K715-TFE 
K816-TFR 
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C CREATE THE NEEDED ELEMENTS OF C MATRIX 
Cl 4— MT 
C25—MT 
C31-CTCE 
C33— STCE 
C66— IG11 
C76— IG21 
C51-STCE 
C53-CTCE 
C86— IG31 
C75-RGP*MP 
C84— RGP*MP 

C GET THE A AND B MATRICIES. SINCE K AND C CONTAIN MOSTLY ZEROS, IT IS 

C BENEFICIAL TO MULTIPLY ONLY THE NECESSARY ELEMENTS INSTEAD OF THE 

C ENTIRE MATRICES 
DO 5 J-1,16 
A(J,1)-Z(J,9)-Z(J,15) 

A ( J , 2) -Z ( J , 10) ~Z ( J , 16) 

A(J,3)-Z(J,ll)-Z(J,14) 

A(J,4)-Z(J,12)-Z(J,1S) 

A(J,5)-Z(J, 13) - Z(J, 16) 

A( J , 6) -Z ( J , 14) *K146+Z ( J , 16) *K166 
A(J,7)-Z(J, 15) 

A(J,8)-Z(J,14)*K148+Z(J,16)*K168 
A(J, 11)-Z(J,6)*K611 

A(J,14)-Z(J,3)*K314+Z(J,5)*K514+Z(J,6)*K614 
A(J,15)-Z(J,4) *K415+Z (J, 7) *K715 
A(J , 16) -Z (J , 3) *K3 16+Z ( J , 5) *K5 16+Z ( J , 8) *K816 
B ( J , 1) -Z ( J , 3) *C31+Z ( J , 5) *C51-Z ( J , 6) 

B(J,2)-Z(J,4)-Z(J,7) 

B (J, 3) -Z (J , 3) *C33+Z (J , 5) *C53“Z ( J,8) 

B ( J , 4) “Z ( J , 1 ) *C 14+Z ( J , 8) *C84 
B ( J , 5) -Z ( J , 2) *C25+Z ( J , 7) *C75 
B ( J , 6) -Z ( J . 6) *C66+Z (J , 7) *C76+Z (J , 8) *C86 
5 B(J,7)— ZCJ.14) 

C SET UP THE INTEGRATION PROCEDURE 
T-TI 

TOL-.0005 

IDO-1 

C SET UP THE PARAM ARRAY REQUIRED BY IMSL ROUTINE IVPRK 
DO 6 J-1,50 

6 PARAM (J) -0.0 
PARAM (4) -3000.0 
PARAM (10) -1.0 

C GET INITIAL VALUES OF INPUT VECTOR U FOR PRINTING AT TIME-0 
CALL UCALC(T,U,AB1,X) 

C SET I.C.'S OF INTEGRAL CONTROL PARAMETER I 

11-0.0 

12-0.0 

13-0.0 
IPT-IPRT 

C SET LIMITS ON INTEGRAL AND ANGLE RATE COMMANDS 
I1MAX-WCI1/KI1 
I2MAX-MCI2/KI2 
I3MAX-WCI3/KI3 
TE1MAX-WCT1/KP1 
TE2MAX-HCT2/KP2 
TE3MAX-WCT3/KP3 

C BEGIN ITERATION LOOP FOR INTEGRATION 

7 TEND-T+DELT 



IH1-I1 

IM2-I2 

IM3-I3 

THETAl “CTCE* (TB 1 +X ( 1 1 ) +X ( 1 4) ) -STCE* (X ( 1 0) +X ( 1 3) +X ( 1 6) ) 

THETA2-X (9) +X ( 1 2 ) +X (1 5 ) 

THETA3-CTCE* (X ( 1 0) +X ( 1 3) +X ( 16) ) +STCE* (TB 1+X (1 1) +X ( 1 4) ) 
IF(T.EQ.TO) GOTO 9 

C IMPOSE LIMIT ON RATE COMMAND DUE TO ANGLE ERROR 
TERR 1 “TC 1 -THETAl 
TERR2-TC2-THETA2 
TERR3-TC3-THETA3 
ATERR1-ABS (TERRI) 

ATERR2“ABS (TERR2) 

ATERR3-ABS (TERR3) 

IF(ATERRl.GT.TElMAX) TERR1“TE1MAX*SIGN(1 . , TERRI) 

IF (ATERR2 . GT . TE2MAX) TERR2“TE2MAX*SIGN (1 . ,TERR2) 

IF (ATERR3 . GT . TE3MAX) TERR3“TE3MAX*SIGN(1 . ,TERR3) 

C IMPOSE LIMIT ON RATE COMMAND DUE TO INTEGRAL CONTROL 
I 1 “ IM 1 +DELT*TERR 1 
12“ IM2+DELT*TERR2 
I3“IM3+DELT*TERR3 
AI1“ABS(I1) 

AI2“ABS(I2) 

AI3“ABS(I3) 

IF(AI1.GT. UMAX) I1-I1MAX*SIGN(1. ,11) 

IF(AI2.GT. I2KAX) I2“I2MAX*SIGN(1. ,12) 

IF(AI3.GT. I3MAX) I3-I3MAX*SIGN (1 . , 13) 

9 CONTINUE 

C NOW FIND THE CONTROL TORQUES 

TP1"KR1* (KI 1*1 1+KP1*TERR1-X (6) ) 

TP2-KR2* (KI2*I2+KP2*TERR2-X (7) ) 

TP3“KR3* (KI3*I3+KP3*TERR3-X (8) ) 

C CONVERT CONTROL TORQUES TO GIMBAL AXIS TORQUES TA,TE,AND TR 
U ( 1 ) -CTCE*TP 1-X (16) *TP2+STCE*TP3 
U(2)-X(16)*CTCE*TP1+TP2+X(16)*STCE*TP3 
U (3) — (STCE+X (15) *CTCE) *TP1+ (CTCE-X (15) *STCE) *TP3 
C IMPOSE TORQUE LIMIT OF TQMAX FT-LBS 
AU1-ABS(U(1)) 

IF(AUl.GT.TQMAX) U(1)-TQMAX*SIGN(1. ,U(1)) 

AU2-ABS(U(2)) 

IF (AU2.GT. TQMAX) U(2) -TQMAX*SIGN (1 . ,U(2) ) 

AU3”ABS (U (3) ) 

IF (AU3.GT. TQMAX) U (3) “TQHAX*SIGN (1 . ,U (3) ) 

IF(IPT.EQ.IPRT) THEN 

WRITE (7) T , X , U , AB 1 , THETAl , THETA2 , THETA3 

IPT-0 

ENDIF 

IPT-IPT+l 

C CALL IMSL ROUTINE TO INTEGRATE 

CALL IVPRK (IDO , 16 , FCN , T , TEND, TOL , PARAM , X) 

IF(T.GE.TF) STOP 

C SET LIMIT ON INTEGRATION STEPS TO AVOID CHANCE OF AN INFINITE LOOP 
IF(T.GE.3600) STOP 
IF(T.LT.O.O) STOP 
GOTO 7 
END 

c 

C SUBROUTINE TO CALCULATE THE R.H.S. OF THE DIFFERENTIAL EQUATIONS 
SUBROUTINE FCN(N,T,X,DX) 

IMPLICIT REAL (K) 
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REAL X (16) ,DX(16) ,U(7) ,A(16, 16) ,Z(16,16) ,B(16,7) ,AX(16) ,BU(16) , 
UtT.MP 

COMMON /SUB1/ A,B,Z , WT,MT,RG,RGP, WP ,U, AB1 ,MP 
C CALCULATE THE TIME-DEPENDENT INPUTS 
CALL UCALC (T,U , AB1 , X) 

C CALCULATE TIME DEPENDENT ELEMENTS IN K AND A MATRICES 
K110-WT+MT*AB1 
K29— WT-MT*AB1 

K7 1 2— RG* (WT+MT*AB 1 ) +RGP* (WP+MP*AB 1) 

K813-K712 

K79-RGP* (WP+MP*AB1) 

K810-K79 
DO 1 J-1,16 

A ( J , 9) -2 ( J , 2) *K29+Z ( J , 7) *K79 
A(J,10)-Z(J,1)*K110+Z(J,8)*K810 
A(J f 12)°Z(J,7)*K712 

1 a(J,13)-Z(J,8)*K813 

C MULTIPLY A*X AND B*U USING IMSL ROUTINE MURRV 
CALL MURRV(16, 16, A, 16,16.x, 1,16. AX) 

CALL MURRV (16,7, B, 16 , 7, U, 1 , 16 ,BU) 

DO 2 J-1,16 

2 DX(J)-AX(J)+BU(j) 

RETURN 

END 

C 

C SUBROUTINE TO CALCULATE TIME DEPENDENT INPUTS 
SUBROUTINE UCALC (T, U , AB 1 ,X) 

IMPLICIT REAL (K) 

COMMON /SUB2/ AB10, AB1 1 , AB20, AB21 ,AB30,AB31 .OMBO.OKBl .OMGAB.TBIO, 
&TB1 

COMMON /TIMES/ TAB11 , TAB12 , TAB21 , TAB22 , TAB31 , TAB32, TOMB1 , TOMB2 
REAL U (7) , X ( 1 6) 

C BALLOON ACCELERATIONS ARE ZERO EXCEPT DURING TIME INTERVAL 
C DEFINED BY TAB1 AND TAB2 
AB1"0.0 
TB1-0.0 
DO 2 1-4,7 
U(l)-0.0 
2 CONTINUE 

IF(T.GE.TAB11.AND.T.LE.TAB12) AB1“AB10+AB11*T 
IF (T.GE. TAB2 1 . AND. T . LE . TAB22) U(4)-AB20+AB21*T 
IF (T.GE.TAB31 . AND. T. LE. TAB32) U(5) -AB30+AB31*T 
C ANGULAR MOTION OF BALLOON 

IF (T. LT. TOMB1 .OR. T.GT. TOMB2) GOTO 5 
U (6) -OMBO+OMB1*T 
U (7) -OMBO*T+ . 5*T*T*OMB 1+OMGAB 
TBl-.5*0MB0*T**2+T**3*OMBl/6.0+0MGAB*T+TB10 
5 RETURN 
END 



DATA FILE 


&INPUT 

TI-O.O, 

TF-15., 

TC1-0.003, 

TC2-0.003, 

TC3-0.003, 

WCT1-. 000645, 

WCT2-. 000645, 

WCT3-. 000645, 

WCI1-. 000645, 

WC12-. 000645, 

WCI3-. 000645, 

TQHAX-10. , 

DELT-0.05, 

IPRT-2, 

TCE-0.500, 

WC-6.28, 

PM1-65 . 0.PM2-65 . 0.PK3-65 .0, 

KGC-1.0E02, 

RP-15.0, 

RG-25.0. 

LC-350.0, 

WG-2205.0, 

HP-1765.0, 

X«16*0.0, 

IP11-79.0, 

IP22- 1830.0, 

IP33-1830.0, 

IP12-0. , 1P13-0. , IP21-0. , 1P23-0. , IP31-0. , 1P32-0. , 
1G1 1-3425.0, 

JG22-1713.0, 

IG33-1713.0, 

1G12-0. , 1G13-0. , 1G21-0. , IG23-0. , IG31-0. , IG32-0. . 
AB10-0.0.AB1 1-0.0, 

AB20-0.0.AB21-0.0, 

AB30-0.0.AB31-0.0, 

OHBO-O.O.OHBl-O.O, 

OMGAB-O.O, 

TB10-0.0, 

TAB11-0.0.TAB12-0.0, 

TAB21-0.0.TAB22-0.0, 

TAB31-0.0.TAB32-0.0, 

TFA-11.0, 

TFE-4.5, 

TFR-4.5, 

LEND 
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